非線形回帰分析　-　10

二次曲線のパラメータの推定

切片，a₀，をシフト

\(\Large \displaystyle y_i = a_0 + a_1 x_i + a_2 x_i^2 \)

切片を変えて傾きを近似すると，このようなグラフを得ることができます．

ここで，S，Se，は，

・残差平方和

推定値からの残差

\(\Large \displaystyle Se = \sum_{i=1}^{n} \left( y_i -\hat{a_0} - \hat{a_1} x_i - \hat{a_2} x_i^2 \right)^2 \)

a₀をシフトさせたときの，推定値からの残差

\(\Large \displaystyle S = \sum_{i=1}^{n} \left( y_i -a_0 - \hat{a_1} x_i - \hat{a_2} x_i^2 \right)^2 \)

であり，傾きの値a₀を，δ，だけシフトさせて，固定し，その際のa₁, a₂の推定値をソルバーで推定しました．

dS，を見ていただけるとわかるように，推定値，Seが一番小さく，ほぼ左右対称に増加していることがわかります．

グラフ化すると，

のように，二乗+定数できれいに近似できます．

二次曲線の近似においてもきれいに近似でき，

\(\Large \displaystyle y = 3.2476 + 1.1736 \ \delta^2 \)

ここで，分散値は，

・分散

\(\Large \displaystyle Ve = \frac{1}{n-3} \sum_{i=1}^{n} \left(y_i -\hat{a_0} - \hat{a_1} x_i \right)^2 = \frac{Se}{n-2} = \frac{3.2476}{6-3} = 1.0826 \)

であり（a₀，a₁，の二つのパラメータが2つあるので，自由度は，n-3），

\(\Large \displaystyle 1.1736 \ \delta^2 = 1.0826 \)

となるδがSEとなるので，

\(\Large \displaystyle \delta^2 = \frac{1.0826}{1.1736} = 0.9224 \)

\(\Large \displaystyle SE_{a_0} = \sqrt{\delta^2} = 0.96042 \)

と推定できます．

・Rによる推定

Rでの近似を行ってみると，

プログラムは，

xx <- c(0,2,3,4,6,9)
yy <- c(5,10,15,25,50,102)
plot(xx,yy)
fm<-nls(yy~a0+a1*xx+a2*xx^2,start=c(a0=5,a1=1,a2=1),trace=TRUE)
summary(fm)

で，結果は，

とドンピシャ一致することがわかります．

Kyplotにおいても，

と同じ結果となります（ということはこの二つのソフトは同じアルゴリズムで標準誤差を計算？）

次は，a₁，をシフトしてみましょう．